Quantum Circa- Rhythms 
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A class of solutions of nonlinear von Neumann equations exhibits rhythmic properties analogous to 
those found for nonlinear kinetic equations employed in phenomenological description of biochemical 
oscillations. As opposed to kinetic equations, typically solved numerically and involving Hill-type 
nonlinearities, our exact solutions are derived by means of soliton techniques and involve a simple 
quadratic nonlinearity. 
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I. INTRODUCTION 

Perhaps the greatest achievement of a pre-quantum 
phenomenology in chemistry was the Mendeleev periodic 
table of elements. The table, in particular, allowed to 
predict new elements with correctly anticipated proper- 
ties. But to really understand the structure of the table 
one had to invent a simply-looking equation 



= Hty) 



(1) 



where \ip) is a vector from a Hilbert space and H a self- 
adjoint operator. 

The Schrodinger equation Q with appropriately cho- 
sen H allowed to understand the origin and types of 
chemical bonds, that is, the deep structures behind the 
Mendeleev table. Although the identification of H for 
a concrete molecule was of crucial importance, there 
were certain features of the abstract equation that al- 
lowed to determine general properties of those systems 
whose dynamics could be governed by (I). For exam- 
ple, a self-adjoint H was related to conservative and iso- 
lated sytems. The conserved quantity is the average 
energy E = (ip\H\ip). Conservative systems that are 
non-isolated (i.e. non-trivially interact with the outside 
world) were described by another simply-looking equa- 
tion 



ip=[H,p] 



(2) 



where p is a trace-class positive operator acting in the 
same Hilbert space and H is the same self-adjoint oper- 
ator. Eq. @ is known as the von Neumann equation. 

is equivalent to Q if p projects on the direction of 
\ip). The invariant (ip\H\ip) is replaced by E = Tr(Hp). 
If \j) is a normalized vector, then the diagonal element 
Pjj = (j\p\j) = TrPjp = pj is a probability associated 
with the proposition represented by Pj = The 
probability is time dependent if [H,Pj] 0. 

Thinking of the present-day chemistry or biochemistry 
one is struck by the ubiquity of nonlinear kinetic equa- 
tions 0, 0- Physically, all kinetic equations describe 
dynamics of probabilities and in this sense are close to 
Typical kinetic equations describe open systems 



and the associated probabilities are related with com- 
muting propositions. Nonlinearity occurs whenever there 
are feedbacks or cooperative phenomena, and details of 
evolution depend on the number of substrates and phe- 
nomenology used to model right-hand-sides of the equa- 
tions. 

Similarly to the catastrophe theory, phenomenological 
modelling is to some extent an art 3j. Having experi- 
ence with nonlinear evolutions one is capable of choosing 
the right-hand-sides of the equations in such a way that 
certain qualitative elements of the dynamics will have to 
occur, and then it remains to play with parameters. Not 
trying to diminish the role of phenomenological modelling 
we have to admit that phenomenology is a two-edged 
sword. On the one hand, it links qualitative analysis 
with quantitative predictions, and often reveals aspects 
that are initially unexpected. On the other hand, how- 
ever, shapes of the equations are dictated by our needs. 
This should be contrasted with Q and (J2J which can be 
derived from first principles. 

A general property of systems described by or J5J is 
the exact linearity of the evolutions they generate. Lin- 
ear evolutions practically do not occur in biological and 
chemical systems, so linearity is in this context a draw- 
back. Sometimes a linear approximation is justifed if 
some paramters are small, but the generic case is nonlin- 
ear. In Ref. |4| we have argued that a candidate equation 
describing a conservative quantum system in a feedback 
with some environment is the nonlinear generalization of 
the von Neumann equation 



ip=[H,f(p)] 



(3) 



where f(p) is an operator function satisfying the "no- 
feedback, no-nonlinearity" property, formally meaning 
that f(p) — p if p is a projector on some direction \ip) 
in a Hilbert space. Eq. is integrable in the sense 
of soliton theory for any /. The time-invariant is now 
E = Tr [Hf(py\. A soliton technique of solving J3J), 
based on Darboux transformations, was introduced in 
[5j, further developed in |fj, and recently generalized in 
y\. Therefore, as opposed to the kinetic equations that 
typically have to be solved numerically (for exceptions 
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cf. H,ilEi)' we can work with exact analytic solutions. 
This allows us to have access to the time scales and pre- 
cision numerically unavailable. 

The set of solutions of contains all the solutions 
of form p — \ip)(ip\ where \ifj) satisfies A new class 
of "self-switching" solutions occuring only for nonlinear 
/ was discovered in [j| and applied in £| in the context 
of morphogenesis. These solutions are characterized by 
three regimes involving a single switching between two 
asymptotically linear evolutions. 

Still, looking for the most striking manifestations of 
biochemical nonlinear evolutions one arrives at phenom- 
ena that are rhythmic. Periodic switching, with periods 
ranging from fraction of a second to years, is encountered 
at all levels of biological organization, and its successful 
modelling in terms of nonlinear kinetic equations is one 
of the great achievements of computational biology (for 
a review cf. ^2). 

The kinetic equations considered in phcnomcnological 
models of biochemical oscillations involve nonlinearities 
of a mixed linear-Hill type. The dynamical systems are 
typically finite dimensional, the number of degrees of 
freedom corresponding to the number of different sub- 
strates one takes into account. The use of Hill func- 
tions of the form Xy,n{x) — x n /{y n +x n ) is motivated by 
phenomenological arguments and their ability to contin- 
uously interpolate between and 1. The latter feature 
allows for a dynamical mechanism of switching, a prop- 
erty whose origin can be traced back to bifurcations of 
nonlinear maps and catastrophe theory. 

Of particular interest are theoretical studies of circa- 
rhythms. Let us recall that circa-rhythms are "classes 
of rhythms that are capable of free-running in constant 
conditions with a period approximating that of the en- 
vironmental cycle to which they are normally synchro- 
nized" 01 • The examples of the rhythms are the cir- 
cadian (24 hours), the circalunar (28 hours), or the cir- 
cannual (365.25 days) ones. What makes circa-rhythms 
interesting from a formal point of view is the interplay 
between the phases where there exists an external forcing 
(light-darkness periods, say), and the rhythmicity that 
sustains even after the external driving is switched off (as 
in experiments in constant darkness). Formally, the two 
cases may correspond to nonlinearities involving (or not) 
an explicit time dependence of some parameters. 

To give an example, the model of circadian rhythms in 
Neurospora [T^l ] involves four variables and two types of 
Hill functions. The rhythmicity occurs in the model even 
without an external driving, but the case with explicit 
time dependence of coefficients is treated as well. Similar 
constructions are given in 0, 0] . 

In the present paper we address the problem of circa- 
rhythmicity from the perspective of Eq. Q. The main 
idea is to consider a composite Hamiltonian system that, 
as a whole, is conservative. There are two sets of degrees 
of freedom corresponding to its two subsystems. When 
we look at certain averages associated with one of the 
subsystems we find that their time evolution can be de- 



termined by an effective density matrix whose dynamics 
is given again by an equation of the form ((2J), but now 
with explicitly time-dependent coefficients in the nonlin- 
ear function f(p). Therefore, the interaction between 
the two subsystems effectively turns one of them into a 
pacemaker that drives the other subsystem and generates 
its evolution with characteristic bursting patterns exibit- 
ing "day" and "night" phases. As opposed to the phe- 
nomenological models found in the literature we do not 
"force" switching by putting Hill functions by hand. We 
take a simple quadratic nonlinearity satisfying the "no 
feedback, no nonlinearity" condition. The motivation for 
this concrete choice comes from its simplicity and expe- 
rience showing that the nonlinear behavior found here is 
generic and analogous to the one found for more com- 
plicated /'s (cf. the example of f(p) — p q — 2p q ~ 1 with 
arbitrary real q analyzed in Q). The solutions we obtain 
are also interesting in themseleves as first examples of 
solutions of J2J exhibiting more than one self-swiching. 

The paper is organized as follows. We begin with an 
example of a nonlinear von Neumann equation © for 
a composite system and time independent nonlinearity. 
We then find its particular solution that is reducible to 
a new von Neumann equation for a subsystem but with 
explicitly time-dependent nonlinearity. We then general- 
ize the example to the class of time-dependent von Neu- 
mann equations and show how to construct their solu- 
tions. Next we explicitly solve an example and illustrate 
its rhythmic properties. 

II. PRELIMINARIES: DYNAMICS IN A 
SUBSYSTEM 

As a simple example consider 

ig = [uj{b\b 1 -b\b 2 ) + {l + X)a)a,Q\-[XJaJ{ 6 W) 
J 

where [frfcjbj,] = 1 = [a, a^] and 

X = h + b\ + b 2 + 4 (5) 
For g 2 = g we find 

ig = [ujJ + a) a, g] (6) 

describing three independent degrees of freedom: har- 
monic oscillation with unit frequency combined with ro- 
tation in a plane with frequency u>. For g 2 ^ g rotation 
and oscillation get nonlinear ly coupled. 

Let us now eliminate the rotation by switching to a 
rotating reference frame 

p = e lulJt ge- iuM (7) 

where the density matrix satisfies 

ip = [(l + X{t))<Ja,p]-[X{t)a)a,f{p)]. (8) 
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Denoting Y = —i(b\ + b\ — b\ — 62) we find 

X(t) = e iuJt I e -' wJi = X cos ut + Y smut 
Y(t) = e lujJt Ye~ lujJt =Ycosojt- Xsiixut 

and, for any t and t', 

[X(t),X(t')] = [Y(t), Y(t')\ = [X(t),Y(t')\ = 0. (9) 

Since X, Y, and a^a commute we can introduce their 
joint eigenvectors 



X\x,y,n) 


= x\x,y,n) 


(10) 


Y\x,y,n) 


= y\x,y,n) 


(11) 


a)a\x, y, n) 


= n\x,y,n) 


(12) 


X(t)\x,y,n) 


= x(t)\x,y,n) 


(13) 


Y(t)\x,y,n) 


= y(t)\x,y,n) 


(14) 


x{t) 


= x cos u)t + y sin tut 


(15) 


y(t) 


= ycosut — x sin cut. 


(16) 



Now take any time-independent normalized vector 

|V>) = J dxdy^(x,y)\x,y) 

and make the Ansatz pit) = \ip)(ip\ ® wit) where w is 
a density matrix acting only on the oscillator degress of 
freedom. Denoting g(p) = p — f(p) = <8> g{w) we 

rewrite © as 

\ip)(i/j\ <g> iw = \ip)(ip\ <S> [a f a, w] (17) 
+ X(t)\ip)(i/j\ <g> a)ag{w) - \ip)(ip\X(t) <g> g(w)a f a 

Taking matrix elements of both sides of l|17J) between 
arbitrary (x,y\ and \x',y') we obtain 

ip{x,y)tl)(x',y')iw = ip(x, y)tp(x' , y')[a)a,w] (18) 
+ tyiz, y)ip(x', y') (x{t)a)ag{w) - x'(t)g{w)a^aj . 

The Ansatz is internally consistent only if i()(x, y) = 
\J 5{x — xo)S(y — 2/0 )■ The diagonal x = x' = xq, y = 
y' = 2/o leads to 

iw = [a/a, w] — xo(t) [a'a, g(w)], (19) 

%o(t) = Xq cosut + yo sinwi, which is an equation of the 
form @ but with / replaced by 

f(w) =w- x (t)g(w) = (1 - x (t))w + x (t)f(w). (20) 

If / satisfies the "no feedback, no nonlinearity" condition 
f(p) = p for p 2 = p, the same holds for the effective /. 

III. QUANTUM PACEMAKER AND 
FEEDBACK 

The pacemaker is an oscillatory system whose state can 
be modified by an external entraining agent, a zeitgeber 



[l2j . Our main interest in this paper is in the free-running 
oscillation of the pacemaker Xo(t), and the circa-rhythms 
it induces in the absence of the zeitgeber. 

The next level is the coupling between the pacemaker, 
which defines the clock "mechanism" , and the observable- 
level circa-rhythms which define the "hands" of the clock. 
We assume a nonlinear feedback between the hands and 
the mechanism, but the dynamics is nondissipative in the 
sense that the energy of the hands averaged over a single 
cycle of the oscillation is constant. 

Motivated by the analysis of the previous section we 
concentrate on the following class of equations 

iw(t) = f 1 (t)[H,w(t)]+f 2 (t)[H,f(w(t))] (21) 

which define the state of the hands of the clock. The 
"no-feedback, no-nonlinearity" condition reads w(t) — 
fi(t)w(t) + f 2 (t)f(w(t)) whenever w{t) 2 = w(t). Now 
assume we know a solution wo{t) of 0- Then 

w{t) = e-^/o/iW^^/^^^etf/o'/iW^) 

is a solution of (|21|l as can be verified by a direct cal- 
culation. The whole problem of solving (|21(l reduces to 
finding a solution of J2J| , which can be performed by soli- 
ton techniques 

Let us consider a simple but generic example where 
f(w) = (1 — s)w + sw 2 and 

fx (t) = 1 + e cos uit, f 2 {t) — —ecosuit (23) 

The parameter s allows us to compare situations where 
the driven dynamics is linear (s = 0) and purely nonlin- 
ear (s = 1), and for the two cases investigate the role of 
the es. Varying s we can also investigate stability prop- 
erties of the rhythms under fluctuations of the feedback. 

IV. HANDS OF THE CLOCK 

The hands of the clock are described by the Hamil- 
tonian H = a^a = X)^°=o n \ n ) ( n \ °f a harmonic oscil- 
lator type. The frequency of the oscillator is equal to 
unity, meaning that this is a reference frequency for all 
the other frequencies found in the system. Modelling the 
hands by a harmonic oscillator is quite natural for ob- 
vious reasons. The quantum oscillator has in addition 
the appealing property of being delocalized. The oscilla- 
tions occur at the level of probabilities in position space 
and what we call "hands" is a kind of a center-of-mass 
quantity. This type of oscillation is what one expects 
in a biological system, since classical oscillators such as 
pendula or springs are too "mechanical" to be of any rel- 
evance. Intuitively, in our model, the hands move at time 
t to the region of space where the concentration of prob- 
ability is the biggest. The probability density may be 
regarded as a measure of state of a certain spatially ex- 
tended object. In the absence of a feedback between the 
hands and the pacemaker, the hands oscillate harmoni- 
cally with their own internal frequency. We shall see later 
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that the nonlinear coupling may practically suppress the 
internal oscillations of the hands in certain intervals of 
time. What will remain are the sudden bursts occuring, 
roughly, with the period of the pacemaker. 

Quantum harmonic oscillator is an infinite dimensional 
dynamical system and, hence, a solution of von Neu- 
mann equations may be characterized by an arbitrary 
number of parameters determining the initial state of the 
hands. There exists a simple trick allowing to construct 
infinitely-dimensional solutions on the basis of a single 
finite-dimensional one. The trick exploits equal spacing 
of eigenvalues of H. With our choice of units the eigenval- 
ues are given simply by natural numbers. Let us divide 
them into sets containing N elements: {0, 1, . . . , N — 1}, 
{N, N + 1, . . . , N + N - 1}, and so on. Each such sub- 
set corresponds to a block in the Hamiltonian, and each 
block can be represented by a N x N diagonal matrix of 
the form 

H k = kl + diag(0, l,...,N-l)=kl + H . (24) 

As a consequence, in each block we have to solve the 
same matrix equation since a restriction Wk of w to the 
fc-th subspace satisfies 



iw k = [H k ,f(w k )} = [H ,f(w k )}. 



(25) 



The job can be reduced to finding a sufficiently general 
solution of a N x N problem. In each subspace we can 
take a different initial condition and a different normal- 
ization of trace. The whole infinite-dimensional solution 
will take the form of a direct sum 



W(t) = ®T=QPkWk{t,Pk), 



(26) 



^2k=oPk ~ !• The fe-th part depends on pk in a compli- 
cated way since the function f(w) is not 1-homogeneous, 
i.e. /(Aw) 7^ A/(iy). The inhomogeneity implies that 
change of normalization simultaneously rescales time; the 
normalization of probability implies that a change of pk 
in a fc-th subspace influences all the other subspaces by 
making their dynamics faster or slower. In this sense the 
solution, in spite of its simplicity, is not a simple direct 
sum of independent evolutions. 

In order to illustrate the possible effects we can use the 
solutions derived in £| for the simplest nontrivial case 
involving self-switching, i.e. for N — 3 and quadratic 
nonlinearity. Of course, since N = 3, the Hamiltonian 
H and the solution w possess matrix elements describing 
transitions 1 <-> 2, 2 «-> 3, 1 <-> 3, between three different 
basis vectors. Although the observation is in itself trivial, 
we want to stress here the formal analogy to the so-called 
three-variable models of bursting oscillations in enzyme 
reactions with autocatalytic regulation 0, . 

We select a subspace spanned by three subsequent 
vectors \k), \k + 1), \k + 2). The family of interest is 
parametrized by a € R controlling the "moment" and 
type of switching between bursts. The parameter nat- 
urally occurs at the level of the Darboux transforma- 
tion, where it characterizes an initial condition for the 



solution of the Lax pair. The density matrices Wk(t) = 
n=o w mn\k + m) (k + n\ are completely characterized 
by the fc-independent matrix of time-dependent coeffi- 
cients w mn . The reader may check by a straightforward 
substitution Il8l that the matrix 



Woo Woi W 2 

W W IVn W12 I = 

W20 W21 W 2 2 



1 



5 m at) \ 
m 5+V5 m 27) 

15 + V5 \ C " W ^ 5 J 



with 



m = (2 + 3*-VK)V3 + V5a cfa , Bt 



V3(ef* + a 2 e-T*) 
9e 2 ^ + (1 + 4V5i) a 



O 2ia>o* 



is indeed a normalized (Triu = 1) solution of the von 
Neumann equation 

iw = [H, (1 - s)w + sw% (28) 
H = diag(0, 1, 2). The parameters are ujq = 1 — ^^ s, 

2 

7 = I5+75 S - 

Now let us rescale the trace. We do it in three steps. 
The modified density matrix 



- p i(l-a)irt, 



>(t)e 



-i(l-s)Ht 



is a solution of iib\{t) — [sH, wi(t) 2 }. Therefore 
w 2 (t) = Ae lil -^ AHt w(At)e- lil -^ Am 

is also a solution of iw2(t) = [H, siU2(t) 2 ] and 
W3 (t) =e-* 1 -^ m w 2 {t)e i{l -^ m 

is a solution of 

iw 3 {t) = [H, (1 - s)w 3 (t) + sw 3 (t) 2 }. 



(29) 
(30) 

(31) 
(32) 



Performing these operations on our explicit solution we 
find 



w 3 {t) = 



A 



5 &(t) ( 3 (t) 
&(*) 5 + Ut) I (33) 



with 

Ut) 

Cs(*) 



15 + VE \C 3 (t) Ut) '5 

(2 + 3i- y/5i) V3 + Tj^ ^A-Ki-sXi-A)^ 
V3(e~< At + a 2 e-~< At ) 
9e 27A t + (i + 4 V5 t ) Q 2 2< ( WoA+(1 _ s)(1 _ A) ) t 

3 ( e 2 7 At +Q ,2) 6 



Now J Q fi(x)dx = t + ^ s'muit, J Q f 2 (x)dx = ——. sinuit. 
illy oh 

w(t) 



We finally obtain the solution 

A / 5 m at)\ 



at) at) 5 / 
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with 

(2 + 3z - y/Ei) \/3 + y/5a 
V3(e-^ A ^ sinwt + a 2 e 7A ^ sinw *) 

xe -i(e W oA+6(l- S )(l-A)-e)^t e ii 

9 e - 2 TA^si„^t + (14. 4V5i) a 2 

CW = 3 ( e -2 7 A^si„ wt+a 2) 

xe -2i(ea) A+e(l- s )(l-A)-e)^£fi e 2it 

of 

Mt;(i) = (l + e cos wi) [H, w(t)] 

-ecosojt[H, (1 - s)w(t) + sw(t) 2 } (35) 

Since TtHw(t) = TriJw(O), TrHw(t) 2 = TriJw(O) 2 , 
the internal energy of the system, averaged over one pe- 
riod T of the pacemaker oscillation, is 

1 r t+T 

E = - I dt'TrHf(w(t')) = TrHw(0), (36) 

and does not depend on i. In this sense the subsystem is 
conservative. 

If one does not integrate over the pacemaker period, 
one finds that the internal energy of the hands harmon- 
ically oscillates with the pacemaker frequency ui. What 
is characteristic, however, the hands do not oscillate har- 
monically but behave as if they were accumulating energy 
during the phases of quiescence in order to suddenly re- 
lease it in violent bursts. Moreover, the bursts are differ- 
ent during the two halves of the pacemeaker cycle, and 
thus resemble the day-night differences one finds in real 
organisms. 

All these properties are illustrated on explicit ex- 
amples. We plot the dynamics of the hands q(t) — 
Trqw(t) = J x ' oo dqq(q\w(t)\q) as functions of time for 
different parameters characterizing the nonlinearity, and 
for different initial conditions. 



timing system whose role is to coordinate various com- 
plicated and different in shape rhythms of a biological 
system. 

Splitting of a "clock" into a "pacemaker" and "hands" 
occurs in our example quite naturally. The entire "clock" 
is a conservative system so the rhythmic behavior we find 
is of endogenous origin. However, the system cannot be 
treated as isolated from external environment since the 
solution p = P®w is not a projector and, hence, involves 
nontrivial external correlations (cf. the analysis of this 
point given in Q). 

To make our analysis simple we have played with a 
concrete family of exact solutions of The family is 
characterized by several parameters related to the choice 
of the equation (to, s) and the initial condition for the dy- 
namics (a^o = — e, a). The solutions we have derived re- 
veal an unknown aspect of nonlinear von Neumann equa- 
tions — the possibilty of evolution involving an infinite 
number of switching events. 

In spite of its simplicity, the quadratic nonlinearity 
generates circa-rhythmic behavior of surprising complex- 
ity and variety. Much more complicated types of behav- 
ior would have been found if one had used the freedom 
of higher dimensionality offered by 126fl . Moreover, the 
solution of the factorized form P ® w, with P a projec- 
tor, is very special and occurs only in cases the operators 
X(t) and X{t') commute for all t, t'. Had we replaced 
the Abelian group of rotations in a plane (with generator 
J) by something more complicated, we would not have 
been able to use the Ansatz with factorization. 

To end these remarks let us make it clear that we do 
not claim to have found a better description of nonlin- 
ear oscillations occuring in biochemical systems than the 
phenomenological approach of Goldbeter and his Brus- 
sels school. What we do claim, however, is that nonlinear 
von Neumann equations of the form (jSJ), in addition to 
their appealing general properties, imply a vast range of 
phenomena with analogies in biophysical systems. One 
should not be surprised if @ will one day play a funda- 
mental role in life sciences. 
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FIG. 1: Circa-rhythm of average position Trqw(t) for e = 2, 
io — 0.08, a = 1, s = —1.1, A = 1. The pacemaker sinusoidal 
oscillation is shown as a reference. 




FIG. 2: The same parameters as in Fig. 1 but with s = —10, 
i.e. for stronger nonlinearity. 




FIG. 4: Trqw(t) for e = 2, a; = 0.9, a = e 2 , s = 1, A = 1. 
Sinusoidal oscillation of the pacemaker shown for reference. 




FIG. 6: Similar situation as in Fig. 5 but for e = 2, ui = 50, 
a = e, s = 30, A = 1. 
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FIG. 7: Stability of circa-rhythmicity under changes of nonlin- 
ear coupling. Average position Tr qw(t) for e — 2, ui = 0.08, 
a = 1, A = 1, and —10 < s < 0. For s sufficiently far from 
the bursts do not qualitatively change with fluctuations of s. 
Fine details of the bursts are smeared out by coarse-graining 
of the plot. 



